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The modulational instability (MI) of spatially uniform states in the nonlinear Schrodinger (NLS) 
equation is examined in the presence of higher-order dissipation. The study is motivated by results 
on the effects of three-body recombination in Bose- Einstein condensates (BECs), as well as by the 
important recent work of Segur et al. on the effects of linear damping in NLS settings. We show 
how the presence of even the weakest possible dissipation suppresses the instability on a longer 
time scale. However, on a shorter scale, the instability growth may take place, and a corresponding 
generalization of the MI criterion is developed. The analytical results are corroborated by numer- 
ical simulations. The method is valid for any power-law dissipation form, including the constant 
dissipation as a special case. 
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I. INTRODUCTION 

The nonlinear Schrodinger (NLS) equation is well 
known to be a generic model to describe evolution of 
envelope waves in nonlinear dispersive media. This equa- 
tion applies to electromagnetic fields in optical fibers [1] , 
deep water waves [2] , Langmuir waves [3] and other types 
of perturbations in plasmas, such as electron-acoustic [4] 
waves and waves in dusty plasmas [5] , macroscopic wave 
functions of Bosc-Einstein condensates (BECs) in dilute 
alkali vapors [6] [in the latter case, it is usually called the 
Gross-Pitaevskii (GP) equation], and many other physi- 
cal systems. The NLS typically appears in the Hamilto- 
nian form 
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where u is a complex envelope field and b represents the 
strength of the nonlinearity (for instance, in the BEC 
context, u is the mean-field wave function and b is pro- 
portional to the s-wave scattering length [6]). The cases 
b > and b < correspond, respectively, to the self- 
focusing and self-defocusing NLS equation. 

In many physically relevant cases, dissipation, which 
is ignored in Eq. (1), is not negligible. In particular, the 
fiber loss gives rise to a linear dissipative term in the NLS 
equation for optical fibers [1,7,8]. Dissipative phenomena 
have also been studied recently in the context of two- 
dimensional water waves [9,10] in order to explain exper- 
imental observations in water tanks. In this setting also, 
the dissipation is accounted for by a linear term. In fact, 
it was the viewpoint of [9,10] on the effects of damping 
that chiefly motivated the present work. Furthermore, 
it has been argued that dissipative phenomena are quite 
relevant in the context of BECs [11,12]. However, in the 
latter case, the loss arises from higher-order processes 
(such as three-body recombination, corresponding to a 
quintic dissipation term in the GP equation [11]) or from 



dynamically-induced losses (such as crossing the Fesh- 
bach resonance [12,13]). Thus, it is necessary to study 
models of the NLS/GP type with more complex forms of 
dissipative terms. Such terms were also introduced in re- 
cent studies of collapse in BECs with negative scattering 
length [14] and evolution of matter-wave soliton trains 
[15]- 

Higher-order loss - in the form of a term accounting 
for five-photon absorption - was also found to play a cru- 
cially important role (along with a cubic term that takes 
into regard the stimulated Raman scattering) as a mecha- 
nism arresting spatiotemporal collapse in self-focusing of 
light in a planar silica waveguide with the Kerr nonlin- 
earity [16]. In this case, the high order of the absorption 
is explained by the fact that the ratio of the resonant ab- 
sorption frequency to the carrier frequency of the photons 
is close to 5. 

One of the most important dynamical features of 
the NLS equation is the modulational instability (MI) 
of spatially uniform, alias continuous- wave (CW), solu- 
tions. The onset of MI leads to the formation of solitary 
waves and coherent structures in all of the above settings 
[1-5,8,9,17,18]. MI is a property of special interest in the 
above-mentioned dissipative models. In particular, as it 
has been recently argued [9,10] in the context of water 
waves, dissipation may lead to the complete asymptotic 
suppression of the instability and of the formation of non- 
linear waves. 

The aim of the present work is to examine the MI in 
the context of the NLS equation with nonlinear dissi- 
pation. Motivated by the model for BECs developed in 
Ref. [11], but trying to keep the formulation as general 
as possible, we thus consider the damped NLS equation 
of the form: 



iu t = —u xx + b\u\ 2 u + (ci - ic 2 )\u\ 2k u, 
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where b and c\ are real parameters. The constant C2 is 
real and positive (in order for the corresponding term to 
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describe loss, rather than gain). The value k = brings 
one back to the much studied case of the linear dissipa- 
tion [1,7-9], while k — 2 corresponds to the quintic loss, 
accounting for the above-mentioned three-body recombi- 
nation collision in BECs [11], and k = 5 - to the above- 
mentioned case of the five-photon absorption, which is 
believed to be the mechanism arresting the spatiotempo- 
ral collapse in glass [16]. 

Our presentation is structured as follows. In Section 
II, we consider the MI in the general equation (2). In 
Section III, we focus on three special cases: (i) the (pre- 
viously studied) linear damping case, (ii) c 2 = [the 
Hamiltonian limit of Eq. (2)], and (iii) the case of purely 
dissipative nonlinearity, with c\ = b = 0. In Section 
IV, we complement our analytical findings by numerical 
simulations. Section V concludes the paper. 



II. MODULATIONAL INSTABILITY FOR THE 
DISSIPATIVE NLS MODEL 



An exact CW (spatially uniform) solution to Eq. (2) 



is 
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where /3o is an arbitrary phase shift, 



F(t) = B 2k + 2kc 2 t 



(3) 
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and B is a positive constant. We now look for perturbed 
solutions of the NLS equation in the form 



u(x, t) — u (t) [1 + ew(t) cos(<7x)] 



(5) 



where w(t) and q are the amplitude and constant 
wavenumber of the perturbation, while e — > is an 
infinitesimal perturbation parameter. Substituting (5) 
into (2), keeping terms linear in e, and splitting the 
perturbation amplitude in its real and imaginary parts, 
w = w r + iwi, we obtain the following ordinary differen- 
tial equations (ODEs) 

w. r = q 2 Wi + 2kc 2 (-F(<)) -1 w r , 

Wi = - (q 2 + 26 (^(i)- 1/fe ) + 2fcci (Fit))- 1 ) w r , (6) 

where the overdot stands for d/dt. These equations can 
be combined into the following linear second-order equa- 
tion for w r , 



w r + 2f(t)w r + g(t)w r = 0, 



(7) 
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Equation (7) can be further cast in a more standard 
form, y + h(t)y = 0, where h(t) = g(t) + f 2 (t), and 

y ee w r {t) exp (/o f( s ) ds ) ■ 

It follows from Eq. (9) that, in the case of the self- 
focusing nonlinearity (6, ci < 0), the initial condition 
«v(0) = (which we will adopt below in numerical simu- 
lations) leads to g(0) < 0, provided that the wavenumber 
is chosen from the interval 



2 k \ 



where 



<f <q 2 CI = -ibBt + k Cl Bt 

+ y](bBl + fcciBg*) 2 + 4fc 2 c3flS*. (10) 

The condition g{0) < implies an initial growth of the 
perturbation, similar to the case of MI in the usual NLS 
equation (1) [1,2,18]. It is interesting to note that, some- 
what counter-intuitively, the interval (10) of unstable 
wavenumbers expands (g 2 r increases) with the increase 
of the loss constant c 2 - This result will be verified by 
means of numerical simulations in Sec. IV. 

However, the perturbation growth can only be tran- 
sient in the presence of the loss, which resembles a known 
situation for the MI in a conservative medium where the 
dispersion coefficient changes its sign in the course of 
the evolution (this actually pertains to evolution of an 
optical signal along the propagation distance in an inho- 
mogeneous dispersive medium) [19]. In particular, even 
when g(0) is negative (i.e., when the condition (10) is sat- 
isfied and we are in the -transiently- unstable regime), 
as t increases, git) increases monotonically as well and 
eventually becomes positive. Thus, we may expect that, 
until a critical time to > 0, the function w r monotoni- 
cally grows but then switches from growth to oscillations. 
This prediction also follows from the consideration of the 
asymptotic form of Eqs. (8) and (9) for t — > oo: as the 
dissipation suppresses fit) and the time-dependent part 
of git), in this limit Eq. (7) gives rise to oscillations of 
w r at the frequency q 2 (and, hence, to stable dynamics). 



III. SPECIAL CASES 



In this section we consider, in more detail, particu- 
lar cases of special interest: (i) k=0 (ii) c 2 = 0, and 
(iii) b = Ci = 0. The first case is the well-known linear 
damping limit (given for completeness). The second case 
corresponds to a higher-order version of the conservative 
NLS equation [see Eq. (1)], while the last one is a case 
of a purely dissipative nonlinearity. 
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A. Linear Dissipation 

In the limiting case k — > 0, Eqs. (2)- (7) lead to the fol- 
lowing results: first, the relevant exponentially decaying 
solution to (2) reads 



, d 2 w r 



dw r 



uo(t) = exp(— c 2 i) exp 



— exp(-2c 2 i) - cii + (3 



(11) 



On the other hand, w r now satisfies the much simpler 
equation 



w" + q 2 [q 2 + 26cxp(-2c 2 t)] w r = 0, 



(12) 



which can be identified, e.g., with the equation (15.2.1) 
of [1] (see also [7,8]). An explicit, though complicated, 
solution to this equation can be found using Bessel func- 
tions [8] and the critical time can be calculated explicitly 
(see e.g., [1]). 



B. The conservative model. 

In the case c 2 = 0, an exact CW solution to Eq. (2) is 

u (t) - A Q cxp [-i{bA 2 + Cl Al k )t] , (13) 

where Ag is an arbitrary constant, and the linearized 
equation (7) assumes the form 



w 



+ q 2 [q 2 + 2(bA 2 + k Cl A 2k )] w r = 0, (14) 



thus giving rise to the MI criterion, q 2 < — 2(feA§ + 
fcciAg ). As it follows from here, MI is possible if ei- 
ther nonlinearity is self- focusing, i.e., either b < or 
ci < 0. If b < and the higher-order nonlinear terms 
are also focusing, then the interval of modulationally un- 
stable wavenumbers grows. On the other hand, if the 
higher order terms are defocusing , then the correspond- 
ing interval shrinks and does so faster for higher powers 
(or BEC densities). In particular, the higher-order self- 
defocusing nonlinearity (k > 1 and ci > 0) completely 
suppresses the MI if the CW amplitude is large enough, 
A 2 ^>\b\/(k Cl ). 



C. The model with dissipative nonlinearity 



In the case of b = c\ = 0, Eq. (7) takes a simple form, 



< + 2 [kc 2 /F(t)\ w' r + q*- Ac 2 k 2 / (F(t)) 



Upon setting r = F(t), Eq. (15) reduces to 



w r = 0. 
(15) 
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which is the Bessel equation. Its solutions are 



wv(t) = A\ J\ 



(16) 



(17) 



where A\^. are arbitrary constants, and J\, Y\ are the 
Bessel functions of the first and second kind, respec- 
tively. According to the known asymptotic properties 
of the Bessel functions, this solution clearly shows tran- 
sition to an oscillatory behavior for large r. Notice that 
the case of b = and C\ ^ is explicitly solvable too, 
through hypergcometric and Laguerre functions, but the 
final expressions are rather cumbersome (and hence not 
shown here). 



IV. NUMERICAL RESULTS 

The above results indicate that the linearized equation 
governing the evolution of the modulational perturba- 
tions can be explicitly solved only in special cases. In 
order to analyze the general case, we verified the long- 
time asymptotic predictions following from Eq. (7) by 
numerical simulations of this equation in the case which 
is most relevant to the BECs. As mentioned above, this 
one corresponds to k = 2, while the other parameters in 
Eq. (2) are set, without loss of generality, to b = —1 and 
Ci = (the latter condition is adopted since, typically in 
the applications of interest, the conservative part of the 
quintic term is less significant than the dominant cubic 
term). Lastly, we fix B a = 1 in the CW solution (4), 
which is also tantamount to the general case, through 
a scaling transformation of Eq. (2), provided that c 2 
is kept as a free parameter. Equation (7) was solved us- 
ing a fourth-order Runge-Kutta integrator, with the time 
step dt = 0.001, and initial conditions w r (0) — 0.05 and 
ii> r (0) = 0. 

Figure 1 illustrates the critical wavenumber bounding 
the MI band (of transient instability, if c 2 ^ 0) as a 
function of the dissipation strength c 2 . The solid line 
is the analytical prediction given by Eq. (10), which is 
compared with points representing numerical data. The 
latter were obtained by identifying the wavenumber for 
each value of c 2 , beyond which there was no initial growth 
of the perturbation w r . Good agreement is between the 
theoretical prediction and the numerical observations is 
obvious. 

Figure 2 shows the difference in the evolution of mod- 
ulational perturbations whose wavenumber is chosen in- 
side (left panels) or outside (right panels) the instability 
band, in the absence (top panels) and in the presence 
(middle and bottom panels) of the quintic dissipation. 
The top left panel shows the case of q — 1.4 for c 2 = 0, 
the top right illustrates the case of q = 1.5 for the same 
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C2, while the middle panels demonstrate the evolution of 
w r for |c 2 | = 0.3, for q = 1.4 (left panel) and q = 1.7 
(right panel). The bottom panels show the evolution 
of the norm in the GP equation (right panel) and the 
spatio-temporal evolution of a function representing the 
perturbation, illustrating (for q = 1) the initial growth 
and eventual stabilization. 



V. CONCLUSIONS 

In this work we studied the modulational instability 
(MI) in equations of the NLS type in the presence of 
nonlinear dissipative perturbations, motivated by the re- 
cent developments in the fields of nonlinear optics, wa- 
ter waves and Bose-Einstein condensates (BECs). For 
general power-law nonlinearities with both conservative 
and dissipative parts, we analyzed the spatially uniform 
(CW) solutions and their stability. We found that the 
dissipation generically suppresses the MI, switching the 
growth of the perturbations into an oscillatory behavior. 
We identified some special cases that can be analyzed 
exactly. Analytical predictions based on the asymptotic 
consideration for t — > oo were compared with numeri- 
cal simulations. A noteworthy prediction (confirmed by 
numerical results) is that the band of wavenumbers for 
which the transient MI occurs expands with the increase 
of the strength of the dissipation. 

As the MI is currently amenable to experimental ob- 
servation in BECs [20] , it would be particularly interest- 
ing to conduct experiments with different densities of the 
condensate. With a higher density, three-body recombi- 
nation effects and the corresponding nonlinear dissipa- 
tion are expected to be stronger, affecting the onset of 
the MI as predicted in this work. 
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FIG. 1. The critical wavenumber q C r [see Eq. (10)] of the 
band which gives rise to the initial growth of modulational 
perturbations around the CW solution. The solid line is the 
analytical prediction (10), and the points show the numeri- 
cal results. The parameters are Bo = 1, b = —1, k = 2 and 
ci = 0. 
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FIG. 2. Top panels display the evolution of the modula- 
tional perturbation, w r , in the absence of dissipation (c2 = 0), 
for wavenumbers that do (q = 1.4, left panel) or do not 
(q — 1.5, right panel) belong to the instability band. The 
middle panels display the same, but in the presence of the 
nonlinear quintic dissipation ( | C2 1 = 0.3), for the perturba- 
tion wavenumbers q = 1.4 (left panel) and q = 1.7 (right 
panel) . The bottom panels show a typical result from the full 
simulation of the Gross-Pitaevskii equation for q = 1. The 
left panel shows the space-time evolution of (|u| — |iio|)/|wo|, 
which evolves in a periodic form, whose amplitude depends 
on w. We can see in the colorbar, the initial transient growth, 
followed by the dissipative evolution. The right panel shows 
the best fit (solid line) and numerical (dashed line practi- 
cally coinciding with the solid one) time evolution of the norm 
P = J \u\ 2 dx. The best fit exponent of log(P) vs. log(l + 1.2t) 
was found to be —0.4994 compared with the theoretical ex- 
ponent of —0.5 (according to Eqs. (3)-(4)). 
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